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1.  INTRODUCTION 


A central  problem  in  target  tracking  is  the  determination  of  position  probabilities 
after  an  arbitrary  time  interval  for  an  assumed  distribution  of  target  course  and  speed. 

A general  treatment  of  this  problem  which  requires  only  that  all  contributing  distribu- 
tions are  discrete  is  given  in  reference  1.  One  consequence  of  the  distribution -free 
nature  of  this  treatment,  however,  is  that  the  procedures  which  result  are  excessive 
in  terms  of  machine  time  or  storage  requirements  where  standard  distributions  are 
involved.  In  such  cases,  probabilities  can  be  generated  by  algorithms  based  on  standard 
mathematical  forms  and  the  storage  and  manipulation  of  large  numerical  arrays  is  thus 
circumvented. 

If  mathematical  expressions  for  both  the  initial  position  distribution  and  that  of  the 
motion  are  available,  the  final  position  distribution  is  expressible  in  terms  of  the  param- 
eters of  these  two  distributions.  Thus,  if  the  mathematical  form  of  the  final  distribution 
can  be  established,  the  computation  of  the  required  probabilities  consists  of  the  following 
steps. 

1 . The  computation  of  the  parameters  of  the  final  position  distribution  from  those 
of  the  constituent  distributions  of  initial  position  and  motion. 

2.  The  computation  of  the  final  position  probabilities  from  the  fully-specified  dis- 
tribution function  so  obtained. 


For  the  standard  bivariate  distributions  under  consideration,  the  number  of  param- 
eters necessary  for  their  specification  is  small,  so  the  extent  of  the  computation  should 
compare  very  favorably  with  that  necessary  to  implement  the  more  general  procedures 
derived  in  the  referenced  memorandum.  If  it  turns  out  that  the  stated  computations 
can  be  based  on  simple,  closed-form  expressions,  this  is  self-evident,  but  unfortunately, 
such  relationships  cannot  be  guaranteed.  We  argue,  however,  that  the  lack  of  a closed- 
form  expression  does  not  necessarily  preclude  a satisfactory  computational  procedure, 
and  that  the  inherent  economy  provided  by  constituent  distributions  with  standard  mathe- 
matical forms  is  not  necessarily  eliminated  by  a transcendental  outcome. 

The  specific  problem  to  be  addressed  is  related  to  current  efforts  to  improve  track- 
ing procedures  in  the  fleet.  For  the  purpose  of  this  research  contribution,  the  problem 
may  be  formulated  as  follows . 


Following  observations  of  a target  of  interest,  the  probability  density  function  of  its 
position  (x,  y)  is  closely  approximated  by  the  bivariate  normal  density  with  mean  (px,  p^) 

variances  and  c^a  , and  correlation  coefficient  p . The  direction  of  target  motion, 

cp  , is  uniformly  distributed  and  independent  of  target  speed.  That  is,  if  V is  the  effec- 
tive speed, 


-1- 


for  0 ^ cp<  2n 


Pvcp(V’ v)  = PV(V)Pa>  M 


1 


-55  Pv™ 


Further,  the  probability  density  of  the  effective  speed  over  a time  interval  t may  be 
approximated  by  means  of  a discrete  density  function,  that  is, 

n 

PV(V)  = £ p.  8<V-V  ) 
i=l 


where  V.,  V , . . .,  V are  known  constants.  Required  are  the  probability  density  func- 
12  n 

tions  of  the  cartesian  and  polar  coordinates,  (u,v)  and  (R,  0)  respectively,  of  the 
position  at  the  end  of  the  time  interval.  The  maximum  target  speed  is  such  that  a flat 
earth  approximation  is  acceptable. 

It  should  be  noted  that  the  derivation  of  the  mathematical  forms  of  these  final  posi- 
tion densities  is  but  an  intermediate  result.  The  computation  of  the  position  probabilities 
with  accuracies  adequate  for  particular  applications  by  means  of  the  facilities  available 
may  require  additional  effort  other  than  at  the  programming  level.  While  the  principal 
emphasis  of  this  research  contribution  is  on  obtaining  relationships  applicable  to  the 
problem  stated  above,  the  advantages  of  a solution  in  terms  of  functions  for  which  com- 
putational procedures  already  exist  have  not  been  overlooked.  A solution  of  this  nature 
based  on  the  formal  solution  derived  in  the  main  text  is  presented  in  appendix  B. 
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2.  GENERAL  FORMULATION 


In  this  section  the  general  formulation  that  will  be  used  is  outlined.  The  succeeding 
section  is  devoted  to  the  specific  problem  just  described. 


Let  (x,  y)  be  the  initial  position  of  the  target  of  interest.  The  joint  probability 

density  function  of  this  position,  p (x, y)  , is  assumed  known. 

xy 

Let  p.r  (V,  cp)  be  the  joint  density  of  target  speed  V and  direction  cp , also  assumed 

V cp 

known.  This  joint  density,  which  describes  the  motion,  may  be  conveniently  replaced  by: 
prcp<r>'!’)=TpV<p(f’'p) 


where  r = Vt  is  the  range  traversed  in  the  time  interval  t . 

Finally,  let  (u,v)  and  (R,  0)  be  the  cartesian  and  polar  coordinates  of  the  target 

position  at  the  end  of  this  time  interval,  and  let  p (u,v)  and  p (R,  0)  be  the  corres- 

uv  K B 

ponding  position  densities. 


To  formulate  the  relationship  of  these  densities  to 
proceed  as  follows . 


Pxy(x,y)  and  pr(p(r,  cp) 


we 


The  target  coordinates  (u,v)  after  time  t > 0 , are  related  to  the  initial  coordinates 
(x,y)  by: 


u = x + r cos  cp  ^ 

v = y + r sin  cp  . 

Taking  x and  y as  fixed,  the  Jacobian  of  the  transformation  is: 
r = f(u-x)a  + (v-y)3}^ 

and  the  conditional  density  p^v(u,v  |x,  y)  is  given  by: 

P„Ju,v  lx*  y>  = f <u'x)3  + <v-y)3}’V  [{(u‘x)3  + (y~v)a  A arctan  pj-l  . 

UV  Lffl  L ll  “X  J 

The  unconditional  density  p^u.v)  is  given  by: 


Puv(u’ V)  = n J 'Z  PUV(U’ V I X’ Y)  Pxy(X’ y)dxdy  ‘ 
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Replacing  puv(u,  v |x,  y)  with  the  expression  derived  above  and  using  the  transformation 

(2.1),  with  u and  v considered  as  fixed,  to  change  the  variables  of  integration  from 
x and  y to  r and  cp , we  obtain: 


If  direction  and  speed  are  independent: 

Pr  (r>  CP)  = pr(r)  pm  (cp) 
rep  r cp 

and 


2 tt  qq 

puv(u'v)  = Jq  p («p)  J"  Vr>pxy(u"r  cos  q>»  v_r  sincp)  dr  dq>  • (2.3) 


If,  in  addition,  p (cp)  is  uniform, 
<P 


Puv(u>  V)  = 2rr  4\"Pr(r)pxy(u'r  C0S  ^ ’ v-r  sin  ^ dr  dcp  * 


(2.4) 


Finally,  if  course  and  speed  are  independent  and  the  former  is  uniformly  distributed  and, 
in  addition,  the  distribution  of  speed  is  discrete,  that  is, 


pv<v>-£pis(v’vi>  ■ 


it  is  readily  shown  that: 


n 


p (r)  = Ep,  5<r-r ) . 

i=l 


r.  = V.t 

l l 


When  these  conditions  all  hold, 


Puv(u*  = !02  So  £ pi  6(r ' “ri>  pxy(u  "r  cos  cp , v -r  sin  cp)dr  dcp 


n 


n 


E P.P.(u.v) 

i=l 


(2.5) 
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where 


X 2 TT 

Pj(u,  v)  = ^ £ p^u-r.  cos  cp , v-r4  sin  cp)  dcp 


Clearly  P^(u,v)  is  the  final  position  density  for  the  special  case  of  the  problem 
under  consideration  in  which 


Py(V)  = 6(V-V.)  , 


that  is,  the  speed  is  assigned  the  constant  value  V.  with  certainty.  It  follows  that  if 
Pj(u,v)  can  be  determined  for  i=l,  2, . . n , the  relationship  (2.5)  provides  the  final 
position  density,  p^u.v)  , for  all  cases  in  which  an  adequate  description  of  the  effec- 
tive speed  is  provided  by  a discrete  density  function. 

Where  the  polar  form  of  the  final  position  coordinates  is  preferred,  the  required 
density  pR(;)(R,  9)  is  obtained  from  puv(u,v)  by  means  of  the  transformation 

R = (u3  +v3)^ 

0 = arctan  v/u  . 


Thus, 

pRp(R,  9)  = RPuy(R  cos  0 , Rsin  0) 
n 

= 5Z  R p.P.(R  cos  0 , Rsin0)  . 
i=l  1 1 


(2.6) 


-5- 


3.  NORMALLY  DISTRIBUTED  INITIAL  POSITION 


For  this  case,  the  density  PXy(x>  Y)  is  bivariate  normal  with  parameters  p^  , 

(7x  , cjy  and  p , and  a solution  to  the  problem  is  obtained  by  evaluating  the  integral 


on  the  right  of  the  expression 
1 2tt 

P.(u,v)  = — f p (u-r.  cos  cp , v-r.  sin  cp)  dcp 
1 2tt  "o  xy  l Y l 

and  applying  (2.5). 


(3.1) 


It  is  clear  that  this  task  would  have  been  simplified  if  the  original  coordinate  system 
had  been  translated  and  rotated  to  provide  a new  system  in  which  the  coordinate  variables 
were  uncorrelated.  This  would,  of  course,  result  in  a final  position  density  with  respect 
to  the  transformed  system,  and  a further  transformation  to  effect  the  restoration  of  the 
original  coordinate  system  would  then  be  necessary.  The  approach  we  shall  take  to 
simplify  the  form  of  the  integral  is  essentially  equivalent  to  this. 

We  introduce  the  transformation: 


u'  = (u  -p^)  cos  3 + (v-Py)  sin  3 
v*  = -(u-p^)  sin)3  + (v-Py)  cos 3 

when  3 satisfies  the  equation 

2poa 


tan  2 3 = 


JL2L. 


(a  +cr  ) {n  ~a  ) 
x y x y 


The  inverse  transformation  is: 

u = p + p'  cos  3 - v'  sin  3 
rx 

v = Py  + p’  sin  3 + V cos  3 . 

The  density  function  of  the  transformed  variables,  P.(u\  v’) 
that  of  the  untransformed  variables,  P^(u,v)  , by: 


(3.2) 


(3.3) 


(3.4) 

is  therefore  related  to 
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P.(u\  v')  = P.  { (p  + u’  cos  0 - v'  sin  0),(p  , + iT  sin  0 + v'  cos  0)  } 
i i x y 

1 2tt 

= — f p f(u  + u'  cos  0 - v'  sin  0-r.  cos  cp)  , 

2tr  Jo  xy  1 rv  r \ 


(p^  + u'  sin  0 + v'  cos  0-r.  sin  cp)  } d cp 


Introducing  the  identities: 


cos  (cp  - 0)  = cos  cp  cos  0 + sin  cp  sin0 
sin  (cp  - |3)  = cos  cp  sin  0 + sin  cp  cos  0 


(3.5a) 


and 


cos  cp  = cos  (cp  - 0)  cos  0 - sin  (9  - 0)  sin  0 
sin  co  = cos  (cp  - 0)  sin  0 + sin  (cp  - 0)  cos  0 


(3.5b) 


and  using  the  latter  pair  to  replace  cos  cp  and  sincp  in  the  integrand  of  the  above  expres- 
sion, we  obtain: 

Pi(u’’  V’)  = 2^  r Pxy  fAi(u’’v'’cP)’  V’’  ^dcp 


where 


A.(u',  v\  cp)  = p + {u'-r.  cos(cp-0)  } cos  J3  - (v'-r.  sin(cp-0)}  sin0 


B^u',  v',  cp)  = p^t  + (u'-r.  cos(cp-0)  } sin/3  + (v'-r.  sin(cp-0)  } cos  0 . 

With  0 chosen  in  accordance  with  (3.3),  this  reduces  to: 


P.(u\  V)  = (2t)  (S, 


r p2n  T tU_riCOS< 

uV  4 exp 2r 

L u 


fu’-r.  cos  (cp-0)  } (v’-r.  sin  (cp  -0)  }a 


2S 


^dco 


"3  - 2n-0  r j f(u’-r  coscp)8  (v'-r  sincp)'-  -| 

- ^ <w  I — +—i — )r 


(3.6) 


where  and  Sv  are  given  in  terms  of  the  parameters  of  the  initial  position  density  by: 
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(3.7) 


S ' S 2(l-pa) 
u v 


[-11  j 1 

, 1 3 4(l-p3)Jf| 

T7  V7  T IT7 

T7  "TV3- 

ay  < ax 

°y  x y ' J 

The  integral  appearing  on  the  right  of  this  relationship  cannot  be  expressed  in 
closed  form,  but  may  be  reduced  to  an  infinite  series  of  Modified  Bessel  Functions . The 
method  used  to  achieve  this  reduction  is  given  in  reference  2 and  is  outlined  in  the  ap- 
pendix. It  gives: 


P^u’.V)  = (2n)"1(SuSv)'^exp(K1-i)  £ (-if  .) cos  (Tky'J  (3.8) 

k=0 


where 


= 1 for  k = 0 
= 2 for  k > 0 


and  the  parameters  q'^  , , w'^  , and  v\  are  as  defined  in  appendix  A. 

This  gives  the  final  position  density  in  the  transformed  coordinate  system.  By 
transforming  back  to  the  variables  u and  v , the  coordinate  system  with  reference 
to  which  the  problem  was  originally  stated  is  restored.  The  required  transformation 
is  given  by  (3 . 4),  namely: 

u = p + u'  cos  j3  - v'  sin  |3 

v = p + u'  sin  £ + v'  cosb 


and 


P.(u,v)  = ^{(u-p^)  cos  ft  + (v-p^)  sin/3  ],  { -(u-px)sin/3  + (v-p^)  cos0}J 

= (2tt)  (S^)  * exp(-qJ)  £ (-l^fi^yF.JI^w^cos  (2kv.) 

k=0 


(3.9) 


where  q.(u,v),  w.(u,v)  and  y^(u»v)  » the  coordinate -dependent  parameters  of  this 
density  function  are  obtained  from  the  corresponding  expressions  for  q'^u'jV')  , 
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w'  (u\v’)  and  v\(u\  v’)  , given  in  the  appendix,  by  making  the  replacements  required 
by  the  transformation,  namely: 


u'  = (u-p^)  cos  3 + (v-Py)  sin0 

v'  = - (u  -p  ) sin 0 + (v -p  ) cos  0 . 

rx  ry 


This  permits  the  computation  of  P.(u,  v)  for  i=l,2 n and  hence,  from  (2.5),  the 

computation  of  the  required  final  position  density  Puv(u» v)  . 

If  the  final  position  density  is  required  in  polar  form,  an  additional  transformation 
is  necessary.  For  any  given  values  of  R and  0 , Pj(R  cos  0,  R sin  0)  is  given  by 
(3.9)  with  parameters  q^RcosG,  rsin0),  w^(R  cos  0,  Rsin0)  and  y.(Rcos  9»  R sin  0)  . 

The  required  final  position  density,  p JR,  0)  , is  then  obtained  from  (2.6). 

R u 
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4.  CONCLUSIONS 


It  follows  from  equation  (2.5)  and  the  discussion  following  its  derivation  that  a solu- 
tion to  the  problem  stated  in  the  Introduction  is  dependent  on  a solution  being  obtained  for 
the  special  case  of  this  problem  where  the  effective  speed  can  be  assigned  with  certainty. 

The  problem  of  obtaining  a solution  for  the  special  case  is  apparent  from  (3.1).  If 
adequate  computational  facilities  are  available,  a direct  numerical  evaluation  of  the 
integral  expression  in  which  p (x,  y)  is  a bivariate  normal  density  of  the  most  general 

type,  is,  at  least,  feasible.  In  the  absence  of  such  facilities,  a reduction  to  an  equiva- 
lent, but  more  tractable,  form  is  mandatory.  Consistent  with  this  reasoning,  an  equiva- 
lent form  has  been  determined  and,  in  a formal  sense,  equation  (3.9)  and  the  discussion 
of  its  parameters,  constitute  the  principal  results  of  this  research  contribution. 

This  concluding  section  would  not  be  complete,  however,  without  a brief  examination 
of  these  results  in  the  context  of  our  opening  remarks.  These  remarks  list  two  steps  as 
being  necessary  for  the  computation  of  position  probabilities  consequent  to  the  mathe- 
matical form  of  the  density  being  available. 

No  special  problem  is  evident  in  implementing  the  first  of  these  steps . Instances 
will  occur  where  the  transformed  coordinate  system  is  adequate  or  preferable  for  the 
treatment  of  the  problem  at  hand.  When  this  is  the  case,  the  parameters  of  the  final 
density  are  directly  provided  by  the  simple  algebraic  formulas  listed  in  appendix  A. 

If  results  are  to  be  referred  to  the  original  coordinate  system  or  if  the  range  and  bear- 
ing of  the  final  position  are  required,  only  the  three  coordinate -dependent  parameters 
are  affected  and  the  modifications  involved  are  readily  implemented. 

If  the  form  of  the  expression  obtained  for  the  final  density  were  entirely  novel,  some 
discussion  of  its  properties  would  be  required  before  the  second  of  the  steps  listed  could 
be  carried  out  in  practice.  In  view,  however,  of  the  striking  resemblance  of  our  princi- 
pal result  to  the  mathematical  form  of  the  Nakagami  Density  Function  (reference  3)  which 
has  extensive  application  in  the  theory  of  tropospheric  scatter,  referral  to  the  literature 
on  this  topic  appears  a preferred  alternative.  There  can  be  little  doubt  that  programs 
that  will  provide  much  of  the  information  that  is  required,  with  only  minor  modification, 
already  exist.  It  is  also  apparent  that  tabulations  of  the  Nakagami  Density  Function  may 
be  used  to  determine  the  final  position  probabilities  subject  to  a minor  preliminary  com- 
putation to  determine  the  appropriate  point  of  entry. 
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APPENDIX  A 


INTEGRAL  REDUCTION 


From  (3.6)  and  (3.7), 

? i 


-a  4 

P.(u',v')  = (2tt)  (SS)? 


U V 


-P  L 


(u'-rjCos  cp)a  (v'-^  sincp)a 


0 


dcp 


where 


Let 
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w;(u,,V)  = ri^+^)i 
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v’^u'.v’)  = arctan 


u 


then 


P.(u*,  v')  = (2tt)  (SuSy)  ^ exp  (-q^)  / 2TT  ^exp  {-F^os  2 cp+  w'.  cos  (cp-y’.)  }dcp 

~P 


From  the  theory  of  Bessel  Functions, 

yx>  = ^ J2tt+c  exp(-x  coscp+ikcp)  dcp 
C 

where  C is  any  real  constant,  and 


exp  (-x  cos  cp)  = £ (-1)  ^(x)  exp  (ik  cp) 

k=  -03 


A-l 


Therefore, 


P.(u\V)  = (2n)‘3(SuSv)'*exp(-q’.) 


yp.)  exp  [2ik  cp  + w'.  cosfcp-v’p  ) dcp 

-|3  k=-<* 


-a  i " k 

= (2n)  (SuSy)  ?rexp(K1'i)  £(-1)  yF^ 

k=-» 


^exp  f2ikcp  + w'.  cos  (cp-y\)  1 dcp 

.g  1 1 

. 00 

= (2tt)'1(SuSv)  ‘exp  <-q’)  £(-l)k  yFi> 

k=-«o 

exp  (2ik  y’.)  t~;  ^ Yi  exp(2ik  cp  + w'  cos  cp)  dcp 

1 n -P-Y*  1 


-J.  00 

= (2n)  (SuSv)"*exp(-q'i)  £ (-l)k  yF^  I^-w^)  exp(2ik  y’ ) 

k=-» 

= (2TT)"1(SuSv)^exp(-q*.)  £ (-l)kyF.)  I^w^)  exp(2ik  y^) 

k=-« 

since,  for  k even,  y-x)  = yx)  • Combining  terms  with  indices  k and  -k 

P.(u\ v’)  = (2n)"1(SuSv)'^  exp«.)  £ (-l)k  e^F.)  I <W ) cos  (2k  y’.) 

k=0 


where 


= 1 for  k = 0 
= 2 for  k > 0 . 


This  is  the  required  result  and  the  relationship  (3.8)  is  thus  established. 


APPENDIX  B 


SOLUTION  VIA  THE  NAKAGAMI  DENSITY  FUNCTION 


The  density  function  referred  to  as  the  Nakagami  Density  Function  in  the  concluding 
section  may  be  introduced  in  the  following  way.  Let  x and  y be  normally  distributed 
with  arbitrary  parameters  p'x  , p'  , a'x  » a » and  p’  , and  let 


X = + y3 

cp  = arctan  y/x  ; 


then 


Pv(x)  = r Pv  m (X » q>)  d«p  = X JTZn  Pv  (X cos  cp , x sin  cp)  dtp  . 


2 TT 


X<p 


xy 


This  relationship  may  be  used  to  define  the  Nakagami  Density  Function.  Although  the 
integral  expression  in  this  form  is  unsuitable  for  computation,  we  may  nevertheless  write 

2 jy 

PY(X  |p'x  , Py  , n’x  » CTy  » p')  = vJQ  Pxy(^cos  cp,  xsintp|p’x,  p*y,  a’x,  a’y.  p’)  dtp  . 

We  make  use  of  this  result  by  substituting  variables  defined  in  section  2 and  in 
the  opening  sentence  of  section  3 as  follows: 

2n 

p (r.  lu-p  , v-p  , a , ct  » P)  = r.  f p (r.cos  cp,  r.sincp  |u-p  ,v-p  , a > a » P)dcp 
‘y  i'  rx  ry  x y iJo  xy  i l T 1 rx  ry  x y 


= r.  f2TTp  (u -r.cos  cp,  v-r,  sin  cplp  , p , a a , o)  dcp 
ro  xy  l i Ynx  y x y 


xy 


= 2TTr.  P.(u,v  |px,  py,  ax,  CTy,  p)  dcp 


by  equation  (3. 1).  Then,  from  (2.5), 


PUV(U'Vk-  V °X’  V P)  = £ PiPi<U-V|MX.  My-  Ox.  Oy,  P) 


1 i Pj 
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r.  x i x y x y 
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This  expression  provides  a solution  to  the  problem  stated  in  the  Introduction,  where 


n 


PyW  = E Pi<v-vi> 

i=l 


r.  = V.t, 

l l 

and  p , p , a , a and  p are  the  parameters  of  the  initial  position  density. 
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